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^ ■ Abstract 

Motivated by recent numerical findings [M. Henkel, T. Enss, and M. Pleiinling, J. Phys. A: Math. 



Gen. 39 (2006) L589] we re-examine via Monte Carlo simulations the linear response function of 



'^ \ the two-dimensional Ising model with Glauber dynamics quenched to the critical point. At variance 



with the results of Henkel et al, we detect discrepancies between the actual scaling behavior of the 
response function and the prediction of Local Scale bivariance. Such differences are clearly visible 



K* ' in the impulse autoresponse function, whereas they are drastically reduced in integrated response 

^^ ■ functions. Accordingly, the scaling form predicted on the basis of Local Scale Invariance simply 

provides an accurate fitting form for some quantities but cannot be considered to be exact. 
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I. INTRODUCTION 

The non-equilibrium collective relaxation of pure systems quenched at or below their 
critical points has been recently the subject of a renewed interest in connection with the 
fact that two-time quantities — such as response and correlation functions — display a 
scaling behavior similar to the one observed in glassy systems (l|, Ij, y]. In addition, at the 
critical point, this scaling behavior is characterized by a certain degree of universality 2|, yj] 
which renders it largely independent of the microscopic details of the systems and allows 
its investigation via suitable minimal models, either on the lattice or in the continuum. 
The former can be easily studied via Monte Carlo simulations (or, in a limited number of 
cases, exact solutions are available), whereas powerful field-theoretical methods have been 
developed for the latter. 

In what follows we will be concerned with the behavior of the two-time response function. 
Consider a lattice model characterized by degrees of freedom a^ with x e Z'^, interacting via 
an Hamiltonian HfjcTx}; {hycit)}] = 7io[{o'x}] — Xlx^x/^-xl^), where hy^(t) is a time- dependent 
external field. In the case of the Ising model Ux = ±1) "^oilcTx}] = ~JJ2ixy) ^xC^y {J > 0), 
the sum being over pairs (xy) of nearest neighboring sites. The stochastic dynamics of 
the model, due to the coupling to a thermal bath of temperature T, can be implemented 
according to different rules and yields a time dependence of the value of the local degree of 
freedom a^{t). The linear two-time response function is defined as 

6{a^{t)) 



Rrit, s) 



6hy{s) 
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where r = |x — y|, (...) indicates the average over the possible realizations of the stochastic 
dynamics, and invariance under space translations has been assumed. Causality implies that 
Rr{t, s) vanishes for s > t and therefore we shall assume t > s in what follows. In addition, 
we shall primarily consider the local response function (also referred to as autoresponse 
function) 

R{t,s) = Rr=o{t,s). (2) 

If the system is in thermal equilibrium with the bath, the dynamics is invariant un- 
der time translations and therefore Rr{t,s) actually depends only on r = t — s, i.e., 
Rr{t,s) = -Rr (t;T). Time-translational invariance is naturally broken if, at time t = 0, 
the temperature of the bath is suddenly changed from an initial value T[ — which we shall 



assume to be high enough to yield a disordered state in the system — to a final value 
Tf 7^ Ti. After the quench, the system undergoes a non-equilibrium relaxation towards the 
new equilibrium state at the temperature Tf, which is attained after a time teq- In some 
circumstances, however, t^q = oo and the dynamics retains its non-equilibrium character for- 
ever. Interestingly enough, a robust scaling behavior emerge in this regime as it is revealed 
by two-time quantities such as the linear response function ([1]). In the paradigmatic case of 
ferromagnetic systems, teq = oo — in the thermodynamic limit — if the quench occurs at 
^f < Tc, where Tc is the critical temperature of the system. In such a case, the behavior 
of Rr{t,s) depends on the relation between t and s. We restrict our analysis to the case 
s,T ^ TmicD where Tmicr IS a microscopic time scale set by the dynamics and such that for 
s,T < Tmicr uou- Universal corrections to scaling appear. In the sector s ^ t ':^ Tmicr, referred 
to as short-time separation regime (STSR), time-translational invariance and time- reversal 
symmetry are recovered in local quantities, resulting in an equilibrium-like dynamics: 

i?(t,s) = i?(<=^)(r;Tf). (3) 

This quasi- equilibrium behavior is quite generally expected in systems with slow dynamics 
(see, e.g., Refs. [l|, y, y]). By contrast, for s ^ r^icr and fixed u = t/s > 1, referred to as 
aging regime, one expects the scaling behavior 

Rit,s) = R^^^\t,s) = s-'-''Fniu), (4) 

which is characterized by the scaling function Fr and the exponent a. These being the 
behaviors in two different regimes, it remains to be clarified: (a) how the crossover between 
the forms ([3]) and (jlj) actually occurs in R{t, s) and (b) how Ff>{u) behaves. 

For quenches at the critical temperature these issues can be addressed via scaling argu- 
ments [21,14] and the renormalization-group (RG) approach 51], which yield 

R{t, s) = Ans-'~''u\u - l)-i-Vij(«) (5) 

where the exponent a is related to the usual static and dynamic critical exponents rj and z 

by 

d-2 + r] 
a = , (6) 

n ^ 

9 is the so-called initial-slip exponent [5|], Ar is a non-universal constant (fixed by the 



1 



condition fniu^ 1) = 1) and fniu) is a universal scaling function such that 

/h(1+)^0 and finite. (7) 



fuiu) can be calculated by means of a variety of field-theoretical techniques. As a conse- 
quence of Eq. ([7]), Eq. ([5]) can be cast in the form 

R{t,s) = R'-''^\T;n)R'^^''\u) (8) 

where 

R^'''\T;n) = Bnr-'-'' (9) 

is the equilibrium critical response function with Bji = Ajifji{l~^) and 

i?(°«^)(n) = n%H///j(l+) (10) 

is the non-equilibrium contribution. In the STSR R^°^^ (u) ~ 1 and the response function ([5]) 
reduces to the form (Q, in agreement with what expected from Eq. ([3]). In the aging regime 
^(off) (^^^ -^ I Q^Yif^ ^]-^g invariance under time-translations, typical of equilibrium, is broken. 
Equation (^ shows that in the quench at Tf = T^ the crossover between the quasi-equilibrium 
/^(cq) Q^Yid the aging R^^^^ behavior occurs multiplicatively via the factor R^"^\ 

Unfortunately, no equivalent powerful methods are available for the analytical investi- 
gation of the behavior of the response function after a quench to Tf < Tc. Nonetheless, 



B y, B, 



the exact results for the 



the available numerical simulations of the Ising model 

Glauber-Ising chain 10|, lll| quenched to Tf = Jjj], and the predictions of approximate 



theories of phase-ordering 13|, [ij, ll5|] suggest that the response function for scalar systems 
in the aging regime can be cast in the form of Eq. (^ where a is a scaling exponent whose 
behavior has been investigated in Refs. [7|, [l6|, with 



fn{u) = {u-iy/'^ . (11) 



3. yb 



(See, however, Refs. [17|, ll8| where fn^u) = 1 has been numerically inferred by looking at 
the so-called thermoremanent magnetization — xity [Oj-^]) ^^^ Xk=o(t, [0,s]), see Eqs. f lT8|) 
and (fT9l) below — for the Ising, 3- and 8-state Potts models.) Here Zg is the growth exponent 
which controls the time-evolution of the typical domain size L{t) ~ t^^^o in systems with 



discrete symmetries (e.g., Zg = 2 for the Ising model with Glauber dynamics) [l9f|. Note that, 
at variance with the prediction of the scaling and RG analysis for the quench to Tf = T^, 
/r(^) vanishes for n ^ 1+. Let us stress that, while for Tf = Tc scaling arguments and RG 
approach predict Eq. ([5]) to be obeyed for every t and s, for Tf < Tc this expression holds 
only in the aging regime. 



Due to the lack of general predictions for the scaling function fR{u) of the response 
function it is particularly interesting the proposal by Henkel et al. — referred to as local 



scale invariance (LSI) [20|, l2l|, |22] — stating that the scaling form of the response function 
is constrained by requiring its covariance under a group of local scale transformations which 
generalize the global transformations (6 > 0) x i— i> 6x, (t, s) ^^ (6^t, Vs) underlying the 
scaling behavior of critical dynamics. 

In this paper we investigate the extent to which LSI describes correctly the scaling be- 
havior of the response function R{t,s), focusing on the two-dimensional Ising model with 
Glauber dynamics after a quench from the disordered state to the critical point. In particular 



we revisit the recent findings of Ref. 23 1. 



The presentation is organized as follows: In Sec. [ITlwe summarize the predictions of the 
different versions of LSI, highlighting their qualitative features and discussing the compar- 
ison with the available analytical and numerical results for a variety of different systems 
undergoing non-equilibrium relaxation after a quench at Tf < Tc. In particular we discuss 
those features of LSI that we shall specifically investigate in the numerical analysis. In 
Sec. Illll we present the results of our Monte Carlo simulations for quenches to Tf = Tc. After 
the description of the numerical procedure, we discuss the data for the impulse response 
function R(t, s) in Subsec. IIII Al and for the associated local and global integrated response 
functions xf^t, s) [cf. Eq. (^] and x^f\t,s) [cf. Eq. ([HD] in Subsec. [mB] and [mCl respec- 
tively, comparing them with previous results and with the predictions of LSI. In Sec. IIVI we 
briefly present the comparison between what is currently known about the response function 
R{t, s), after quenches to Tf < Tc, and LSI. We point out that the aging part of R(t, s) is ac- 
tually described by a recently proposed version of LSI. A final summary and our conclusions 
are presented in Sec. |Vl 

II. LOCAL SCALE INVARIANCE: SCALING FORMS FOR THE RESPONSE 
FUNCTION 



For the impulse autoresponse function R{t, s), LSI in its original form 20|, l22l] (which we 
refer to as LSI.O) predicts {u = t/s > 1) 

i?(t,s) = Ajjs"^-°u^(n-l)-^~" (12) 



and 



Rr{t,s) = R{t,s)<^{r{t- s) 



-l/^^ 



(13) 



where the exponents a, (3 and the amphtude Af{ are undetermined whereas $ is defined 
in terms of special functions in Ref. |22(|. These predictions are supposedly quite general 
and they apply to a variety of systems undergoing non-equilibrium relaxation, either due to 
critical (Tf = Tc) or phase ordering (Tf < Tc) dynamics. As in the case of dynamics at the 
critical point (see Eqs. ([5]) and ([7])), the structure of Eq. (fT2l) is such that time-translational 
invariance is recovered in the STSR: R{t = s + r, s ^ r) = A^r"^"" and this stationary 
part crosses over to the non-equilibrium behavior via a multiplicative factor as in Eq. (jHj). 
Note that, according to LSI, such a multiplicative structure is expected in ferromagnetic 
systems both for quenches to Tc and below Tc- The prediction ( IT2l) (with Ar, a and (5 as 



fitting parameters) has been tested in a variety of different systems 



27 



28 



29 



30 




20 



24 



25 



m 



311 , with the conclusion that, while it is seemingly successful in a restricted 



number of instances, it surely fails in the general case. For example, LSI.O is in agreement 



10 



11 



11 



32| quenched to Tf < Tc, 



13l . Il4l . Il5l l obtained for 



with the exact solution of the large-A^ (spherical) model 
but fails to reproduce all the available analytical results 
scalar systems quenched to Tf < Tc- In addition, for the d- dimensional Ising model with 
Glauber dynamics quenched to Tf = T^, a field-theoretical calculation in the dimensional 
expansion d = A — e clearly shows corrections of order e^ to the prediction (TT2l) of LSI.O 3J] 
for d < A. The comparison of LSI.O with Monte Carlo simulations, at and below Tc-, is more 
controversial since agreement was reported in the study of various systems, such as the two- 
and three-dimensional Isin g rn odel 171 . |20| , the XY model 33| and the 3- and 8-state Potts 
model in two dimensions 18|, while deviations from Eq. flT2|) were detected in the two- 



and three-dimensional Ising model after a quench to the critical 
below it in Refs. [Tl, Isj. In particular, the data presented in Ref. 
intermediate response function 



point in Refs. 



and 



29| for the so-called global 



Js/2 „ 



(14) 



and in Ref. 



9| for Rr=o{t,s) display deviations from LSI.O which are actually negligible in 
d = A but become significant in d = 3 and even larger in d = 2, i n q ualitative agreement 
with the expectation based on the field-theoretical analysis of Ref. 34l |. 
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Motivated by these discrepancies, Henkel et al. have reconsidered the original formulation 
of LSI, coming up with a new version 23|, |28| — referred to as LSI.l in the following — 
which predicts, in comparison to LSI.O, a less constrained scaling form 



R{t,s) = ARS- 



l-o^/3+a'-a| 



U 



-1-a' 



(15) 



where a and a' — yet undetermined — are not necessarily equal. Note, however, that the 
large-M behavior of Eq. fllSp . being actually independent of a', is the same as the one of 
LSI.O (see Eq. flT^ ) and therefore possible improvements of LSI.l compared to LSI.O are 
restricted to the STSR m ~ 1. (The space-dependence oi Rrjt, s) within LSI.l still factorizes 



according to Eq. (TT3l) . with R{t,s) given by Eq. (fTSjl [2ll . l23|.) The introduction of the 
additional fitting parameter Aa = a — a' in Eq. flTSjl is expected to reduce discrepancies in 
those cases in which they where observed with Aa = 0, i.e., when comparing with LSI.O. 
A first apparent confirmation of the predictions of LSI.l came from the comparison with 
the numerical studies of (i) the two- and three-dimensional critical Ising model (Glauber 



dynamics) [23|] , (ii) the three-dimensional Ising spin glass |23l . ISSL (iii) the non-equilibrium 
kinetic Ising model in one dimension 26|, (iv) the contact process J23|, and with the analytical 
study of the Fredrickson- Andersen model [23[ . However, subsequent analyses of the latter 



two cases have revealed clear discrepancies with LSI.l |30l. I3ll.l36|. 



In this paper we revise the statement of Ref. 23|] according to which the numerical data 



,(g) 



for Xi {t,s) (see Eq. (HI 



in the two-dimensional Ising-Glauber model — though not in 
agreement with LSI.O [29| — do scale according to LSI.l flTSj) with a = a — d/z, (3 = 
(see after Eq. ([5])) and Aa = 0.187(20). More specifically, we shall focus on the fact that 
the prediction flTSl) of LSI.l with Aa ^ implies a non- stationary behavior of the response 
function in the STSR: 



i?(t = r + s, s > r) 



Ars-^^t"^-"' 



(16) 



which is in qualitative disagreement with what generally expected (see Eq. (|3])) also on the 
basis of the scaling and RG analyses (Eqs. ^ and ([7])). A conclusive test of LSI.l for Tf = T^ 
can therefore be done by looking at the STSR. Accordingly, we performed new extensive 
Monte Carlo simulations in order to investigate the scaling behavior of the response func- 
tion deep in the short-time separation regime, down to r/s ~ 5 ■ 10"^. The numerical data 
presented in Sec. IIIII unambiguously show that R{t, s) becomes time-translational invariant 



in the STSR, in agreement with Eq. ([3]) and in stark contrast with the prediction (fT6|) (with 
Aa 7^ 0) of LSI.l. Nonetheless, the data for X/ (^) ■5) reproduce the findings of Ref. 23| 



showing thereby that in passing from R{t,s) to Xi {'t,s) the space-time integration drasti- 
cally reduces the difference between the actual response function and the prediction of LSI.l. 
We corroborate this statement by studying the scaling behavior of the local intermediate re- 
sponse function 

X?{t,s)= r dt'R{t,t'), (17) 



Js/2 

which involves only one integration over time and which actually deviates from LSI.l less 
than in the case of R{t,s) but more than Xi {^,3). Therefore, in spite of the agreement 
of X/ (^) ■5) with LSI.l, we conclude that the prediction of Local Scale Invariance simply 
provides an accurate fitting form for some quantities but cannot be considered to be exact 
for Tf = Tc. 

III. QUENCH TO T = T^: MONTE CARLO SIMULATIONS 

We consider the Glauber dynamics of the ferromagnetic Ising model on a two-dimensional 
square lattice with 10'^ x 10'^ spins. The system is initially prepared in an infinite-temperature 
configuration and then it is quenched to the critical temperature Tf = Tc ^ 2.269185. 
Temperature is measured in units of the ferromagnetic coupling J and time is expressed in 
Monte Carlo step units (sweeps). Each of the data points we shall present in the following 
is the result of an average over 10^ — 5 ■ 10^ different realizations of the initial condition and 
of the thermal noise. No finite-size effects were detected in the time range accessed during 
the simulations. 

For later convenience, we introduce the (time) integrated response functions in real and 
momentum space 

X{t,[h,t2])= f'dt'R{t,t') (18) 

and 

Xk{t,[ti,t2])= f'dt'Rk{t,t'), (19) 

Jti 

respectively, where Rk{t,s) = Ylr^r{t,s)e''''^'^ is the Fourier transform of Rr{t,s). The 

various response functions we are interested in, namely i?(t,s), X/ (^^s), and X/ (^) ■s) are 

related to xi^^ [^1)^2]) and Xfc(^) [^1)^2]), which are actually computed without applying the 
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external perturbation hx{t) via the algorithm presented in Ref. [37|]. In fact, these response 
functions can be obtained from particular correlation functions of the unperturbed system 

as 



37|: 



2T,x(t, [ti,t2]) = C{tM) - C{tM) + D{t, [t^M]) (20) 

and 

2reXfc(t, [tiM) = Ckit, ti) - Ckit, ta) + Dk{t, [ti,t2]), (21) 

where C(t,t') = (crx(t)crx(i')) is the spin autocorrelation function, (...) being the average 
over initial conditions and thermal histories, and Ck(t,t') = (o'k(t)(T-k(^')) is its momentum 
counterpart. The additional correlation functions D and D in Eqs. f l2D]) and T2Th are re- 
spectively given by D{t, [ti,t2]) = Mt)B^{[ti,t2])) and Dk{t, [h,t2]) = (<Tk(t)fi_k([ti, ta])) 
where By_{t, [^1,^2]) is a function of the transition rates and B^it, [ti, ^2]) is its Fourier trans- 



one has 

"t2 



form. For Glauber dynamics [S 

B^iihM) = r dt' [a^it') - tanh (i/f (t')/T,)] , (22) 

Ju 



where H^{t) is the Weiss field H^{t) = Xl(xy) '^y('^)' ^^^ sum being over sites y nearest 
neighbors of x. Note that the (auto) response x and the correlation C appearing in Eq. (l20|l 
do not depend on x due to space homogeneity. In order to reduce the noise of the resulting 
numerical data, they have been computed at each lattice site and then averaged over the 
whole lattice. 



A. Impulse autoresponse function 

The impulse autoresponse function R can be formally derived from x(t, [^1,^2]) (see 
Eq. ( ITSj) ) as R{t,s) = lim.t^-.sXi'tj [■s,t2])/(^2 — s). Numerically, R is computed as R(t,s) = 
{l/6)x{t, [s,s + 6]), where the integration over time in x has the effect of reducing the noise- 
to-signal ratio with the drawback of introducing a systematic error of order S/s |37|. With 
an appropriate choice of s and 6, this error can be made much smaller than the statistical 
errors and hence neglected. In the range of s considered in our simulations (s > 10^), 6 = 1 
turns out to be a suitable choice. 

The behavior of R{t, s) for t ^ s, i.e., in the aging regime, is quite well established and 



numerical results for x(t, [0, s]) (see, e.g., Refs. [4, l20||), R{t,s) [9|, l39||, and Xj (t, s) [29| 



M)i 



support the validity of the scahng behavior predicted within the scahng and field-theoretical 
approaches (see Eqs. ([5]) and ([7])): R{t ^ s, s) ~ s~^~"(t/s)^ where the actual exponents are 
compatible with the expected values in d = 2 (see Eq. (^): a = rj/z ~ 0.115 {r] = 1/4 and 
z = 2.1667(5) 40|) and 9 = 0.383(3) [3|, |4l|. Accordingly, the comparison of the observed 



behavior of R{t, s) in the aging regime with LSI.l fixes a = a and f3 = 6 in Eq. (1151) . leaving 



Aa as the only undetermined exponent. In Ref. 23 1 



Aa = 0.187(20) (23) 

has been determined as the best parameter value to fit — actually with no visible corrections 
— the numerical data of Xi i'^^s) with the prediction flTSj) of LSI.l. As pointed out in 
Sec. [ni if the scaling behavior of R{t, s) is actually captured by LSI.l with Act ^ then in 
the short-time separation regime t ^ s ^ Tmicr one will have R(t, s) ~ ^-Aa^-i-a+Aa ^g^g 
Eq. (IT6l) ). instead of the generally expected quasi-equilibrium behavior R{t, s) ~ r~^~" (see 
Eq. ([3]) and (Q). In order to test these two qualitatively different predictions we focus on 
the behavior of TcR{t, s) in the STSR, reported in Fig. [H The scaling behavior is expected 
to set in for t,s ':^ Tmicr, i-e., the actual STSR is restricted to Tmicr <C r ^ s. In the 
inset of Fig. [1] we report, for fixed r = 10 ^ Tmicr — 1, the behavior of TcR{t,s) upon 
increasing s, eventually probing the STSR. In particular, according to LSI.l, one expects 
R{s + 10, s) ~ s~^", whereas the actual data seemingly attain a plateau for s > 3 ■ 10^ 
(highlighted by a box in the inset) setting an upper bound to the asymptotic value of Aa: 
Aa < 0.029(8), which definitely excludes the value (!23l) indicated in the inset as a dashed 
line. It seems therefore reasonable to conclude that Aa = for large enough s. Note, 
however, that the behavior of R{s + 10, s) for s < 10^ could be easily fitted by a power law 
with an effective exponent Aa^s — 0.1 which seemingly increases upon decreasing the largest 
value of s included in the fit. If this observation carries over to the case of Xi (ty s) analyzed 



in Ref. (23|], it might explain the quite large Aa (1231) which has been determined there on 
the basis of data with s < 200. Repeating the previous analysis for different values of r we 
conclude that for s > 3 ■ 10^ (highlighted by a box in the inset) and r < 10^ (corresponding 
to t/s < 0.032), R{t + s, s) is actually very close to its plateau value in the STSR, although 
small corrections of order t/s are still present. The data corresponding to these selected 
range of r and s, reported in the main plot in Fig. [H display a quite good collapse onto a 
master curve R{t+s, s) ~ r^^~"'^'^ with Cmc — 0.095(2) when fitted in the region 5 < r < 30 
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FIG. 1: Autoresponse function R{t,s) as a function of r = t — s for different values of s = 
4 ■ 10'^, 5 • 10^, 7.5 • 10^, and 10^. The dashed hne indicates the power-law behavior R{t, s) ~ r~^~" 
with a = 0.115 (see Eq. (HI)). The box indicates the points used to fit the value of the exponent 
flMC (see the main text). In the inset R(t,s) is plotted as a function of s for fixed r = 10. The 
dashed line indicates the behavior oc s~^°' predicted by LSI.l, with Aa given by Eq. ([25]) . The 
box contains the points also present in the main figure. 

(the corresponding points are highlighted by a box in the main plot of Fig. [T]). This value 
is in reasonable agreement both with a = 0.115(2) and with a — Aa, which are expected 
on the basis of Eqs. ([3]) and iQ and of LSI.l with Aa < 0.02, respectively. In comparing 
these values one has to take into account that our Monte Carlo estimate omc — 0.095(2) is 
seemingly affected by a systematic correction which reduces omc compared to its asymptotic 
value for s — i> oo. 

Although the analysis of the scaling behavior of the response function in the STSR does 
not allow us to rule out the validity of LSI.l, at least it definitely indicates that Aa < 0.02. 

A more effective comparison between the numerical data and the different theoretical 
predictions can be done by considering the quantity 



gn{t,s)=T'+''u-'R{t,s 



(24) 
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for arbitrary t > s. According to the general scaling behavior Eq. ([5]), one expects 

gR{t,s)=ARfR{u), (25) 

i.e., gji(t,s) depends only on u = t/s and is such that limj;^., 5f^(t, s) = Ar. In addition, 
within the RG approach, fni^u) satisfies the condition (ITj), whereas within LSI.l fR^u) is 
given by 

/^lsi-i)(m) = m-^"(m - 1)^°, (26) 

(compare Eq. (llSp with Eq. (^) resulting in two different limiting behaviors for m — > 1: 

Br, RG, 

\imgR{t,s)={ (27) 






Aj,(u-1)^", LSI.l, 



where Br = Ar/r^I^) ^ (see after Eq. ([9])). The qualitative difference between RG and 
LSI.l in the STSR translates into the fact that for small values of the abscissa they predict 
the occurrence of a plateau and a power-law behavior with exponent Aa, respectively, when 
5'/j(t, s) is plotted as a function of n — 1. The log-log plot of gRit, s) is presented in Fig. [2] for 
s = 10^, 2 ■ lOf, 7.5 ■ 10^, and 10^, which extend considerably the range of times investigated 



in Refs. 



23 



29[, with s < 200. We report only data with r > 2 because, as discussed in 
Sec. m scaling sets in only for r ^ Tmicr — 1- For a fixed value of s, the statistical uncertainty 
affecting the data points in Fig. [2] increases upon increasing u — 1, i.e., t, due to the fact 
that the value of the response function R{t, s) ^ becomes increasingly small and therefore 
comparable with statistical fluctuations. As a consequence of the typical scaling (jll), this 
effect is expected to be more severe upon increasing s, as it is conflrmed by comparing the 
data points with s = 10^ and s = 10^. The resulting data collapse in Fig. ([2]) is quite good 
(within the errorbars) and conflrms that no signiflcant corrections to scaling are present. 
From the behavior of gR(t, s) at large u — 1 one determines the non-universal constant Ar. 

Ar = 0.071(2), (28) 

(obtained by fltting the data with s = 100 and u > 7) whereas in the opposite limit u ^ 1~^ 
one flnds (from the data with s = 7.5 ■ 10^, 10^ and m — 1 < 1.2 ■ 10"^): 

Br = 0.040(1), (29) 

and therefore one estimates 

/«(1+) = Br/Ar = 0.56(3). (30) 
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FIG. 2: Scaling function gji{t, s) (see Eq. (124p ) as a function of u — 1 (u = t/s > 1) for s = 100, 
2000, 7500, and 10000. The solid line is the prediction of LSI.O whereas the dashed and dash-dotted 
lines are those of LSI.l with Aa = 0.187 and Aa = 0.115, respectively. In all these cases the values 
9 = 0.38, a = 0.115, and Afi = 0.07 have been used. 

1. Comparison with the available field-theoretical results 

Note that gR{t, s)/Ar = JrIu) characterizes the universal scaAing behavior of the response 
function within the Ising universality class with purely dissipative dynamics. In turn, this 
fact can be used in order to determine the associated scaling function in simplified models, 



e.g., field-theoretical (FT) ones 



in Ref. 



421 ] belonging to this universality class. In particular. 



42l | the scaling function for the spatial Fourier transform of Rr(t,s) (Eq. ([T])) was 



determined up to the first order me = A — d>0 where d is the spatial dimensionality of 
the model. This result allows the calculation of fR^u) (see Eqs. (^ and (^^), yielding 



/r'l 



u 



l + 0{e' 



(31) 



which actually does not reproduce the non-trivial behavior observed in Fig. [2] and the actual 
value of /ij(l"^) determined numerically (see Eq. (1501) ). On the other hand, as in the case of 
the scaling function of "^^.Rrit, s) studied analytically in Ref. 3J|, one expects non-trivial 
corrections to /jj (u) at O(e^) which might reproduce at least the qualitative features of 
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the numerical data (such as the fact that gR(t, s) < A^). The calculation of such corrections, 
however, is beyond the scope of the present study. 

2. Comparison with LSI 

The mastercurve in Fig. [2] clearly exhibits a plateau at small u — 1 which excludes the 
power-law behavior predicted in this regime by LSI.l (or, at least, it reduces significantly 
the upper bound to Aa which results from Fig. [1] 43|]). This conclusion is also confirmed 
by comparing the numerical data for gR{t, s) with the prediction fl2Sl) . where fR{u) is given 
by Eq. fl26l) and Ar by the fitted value fl28l) . In Fig. [2] we present this comparison for 
Aa = (solid line), corresponding to LSI.O, Aa = 0.187 (dashed line), corresponding 



to the estimate fl23l) of Ref. 23|, and Aa = 0.115 (dash-dotted line) — which we have 



determined from the best fit to the actual data. Clearly, both LSI.O and LSI.l with the 



value of Aa determined in Ref. [23|] do not capture correctly the scaling behavior of the 
response function already for n ~ 6 and in particular they overestimate and underestimate 
the numerical gR{t,s), respectively. The choice Aa = 0.115, however, results in a good fit 
of the numerical data down to n — 1 ~ 2 ■ 10~^. In spite of that, sizable discrepancies are 
anyhow observed for smaller values of m — 1 (STSR) and no further improvement is actually 
possible by varying Aa. In passing we mention that a similar situation is also encountered 



when comparing the predictions of LSI.l with the scalin g be . 
function within the directed percolation universality class 30 



lavior of the aging response 



3l|. 



The conclusion we can draw so far is that the scaling behavior of the impulse response 
function R{t, s) does not agree with the prediction of LSI.l even though the latter actually 
provides a good fit to the numerical data for t/s — 1 > 2 ■ 10~^. 

B. Integrated local response functions 

In this Subsection we discuss the behavior of the local integrated response function (see 
Eq. (TT71)). which is given by Xi i'^i ^) = x(^? [-5/2, s]) in terms of the quantity ( ITSl) measured 
in the simulations. This integrated response function has also been considered in Refs. [7|, |8| . 

According to the general scaling behavior of R{t, s) (see Eq. ([5])), one expects 

xf{t, s) = ARS-'^u'-'iu - l)"Vx,i(^) (32) 

14 



where the universal scahng function fx (u) is given, in terms of the scahng function fji, by 
/^,i(w) = il-u-'r dx x-^ (l - -) " "" fniu/x) , (33) 

J 1/2 ^ ""^ 

and is such that fx,iiu ^ I) = Kg with kq = {1 — 2^-'^+^)/(l — 9) {kq ~ 0.564 in two 
dimensions). In particular, the prediction of LSI.l is obtained from f\^ ' [u] given in 
Eq. ([26]): 



/x^'^H = ^-^-[^r~ ^'^'^^ - ^, 1 + « - A«, l/u) - 2'~\F,{1 - ^, 1 + a - Aa, l/(2n))] , 

(34) 
where 2-^1 is the hypergeometric function and the case of LSI.O is recovered for Aa = 0. In 
order to compare the different theoretical predictions we consider the function 

g^,,{t,s)^r'^u'~\f{t,s) (35) 

which, according to Eq. ( l32l) . is a function of m = t/ s 

9xAt^ s) = Anf^^iiu) (36) 

such that A^^i = limj>sfl'x,i('^) s) = ArKq. In the opposite limit, instead, 

{B-^i = Bj^/a, RG, 

(37) 
A^(n-l)^"/(a-Aa), LSI.l, 

which is valid under the assumption a, a — Aa > and which reflects the different behaviors 
of /i?(n -^ 1+) within RG and LSI.l (see Eqs. (171) and (l26l) . respectively). In Fig. [3] we 
plot q^ ] (t, s) for different values of s, which are slightly larger than those considered in 



Ref. 



23 



291 but smaller than the larger values which Fig. [2] refers to. For s > 600 and 
M — 1 > 10"^ one observes a rather good collapse onto a single master curve which depends 
only on u, in agreement with the scaling (l32l) . At smaller values of m — 1 the actual scaling 
behavior is accessed only for increasingly larger values of s such that t = t — s > Tmicr, i-e., 
for s > Tmicr/{u — 1)- From the large-M behavior of gx,\{t, s) we determine the non-universal 
constant 

A^i = 0.036(3) (38) 

in marginal agreement with the expected value Af(Kg = 0.040(1) resulting from the estimate 
of Ar given in Eq. (I28l) . Note that while gR{t, s) < Ar, here we find gx,i{t, s) > A^^y. 
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FIG. 3: Scaling function g^^i{t,s) (see Eq. (|35p ) as a function of u — 1 (u = t/s > 1) for s = 50, 
200, 600, and 2000. The solid line is the prediction of LSI.O whereas the dashed and dash-dotted 
lines are those of LSI.l with Aa = 0.187 and Aa = 0.115, respectively. In all these cases the 
values 9 = 0.38, a = 0.115, and Ajj = 0.07 have been used. The horizontal dashed line on the left 
(right) indicates the limiting value Bn/a {Arkq) expected for n ^ 1 (n — > oo) according to the 
RG and scaling behavior (see Eq. (j37p ) where Ajj and Bji have been determined from the behavior 
of gR{t,s) (see Eqs. ([25]) and 



1. Comparison with the available field-theoretical results: 



The universal scaling behavior of X/ (^) s) within the Ising universality class with purely 
dissipative dynamics is characterized by g y.i (t, s) /Ay,,] = fx,\{u) which can also be studied 
within field-theoretical models [3|, |3J, |4^. As a consequence of the fact that /)j (u) = 
1 + 0(e2) = /]f^^-°V) + 0(e2) (see Eq. (EH) and Eq. ([2SD with Aa = 0), one finds, from 



[■(LSI.O), 



(.(LSI.O), 



Eq. dMD, 47V) = fxT^'^'i^) + 0{e^), where /^7"''^(m) is given by Eq. ([MD with Aa = 0. 
The resulting prediction ^_r/^i (m) for gx,i{t, s) is shown in Fig. [3] as LSI.O: The differences 
with the numerical data are less severe than in the case of A^f}^ (n) for gnit, s) in Fig. [2l 



In particular, for m — i> 1+ one has //, [ (l"*") = /. 



(LSI.O) , 
'X,l ' 



-0(e2) = l/a + 0{e^) (see Eq. §^ 



with Aa = and u -^ 1"*") 4J], yielding the plateau in the FT prediction reported in Fig. [3] 
as a dashed line, which overestimates the actual numerical value. It would be interesting to 
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see whether the corrections of O(e^) to /^j (n) account for such a difference. 

2. Comparison with LSI: 

As in the case of the impulse response function, the different quahtative behavior of the 
LSI.l and RG predictions should be displayed for u — 1 ^0 according to Eq. fl37j) . but 
the data presented in Fig. [3] are not actually able to detect such a difference. However, the 
trend of the data is compatible with the RG prediction Urates gx,\{t, s) = Br/ a ^ 0.35 (see 
Eq. (l29l) ). For comparison we plot also Eq. (l36l) with K\ predicted by LSI (see Eq. ( IMl) ). 
An = 0.063 (note that this is slightly smaller than the previous estimate Eq. ( 128|) ) and for 
the different values of Aa already considered in Fig. [2j Aa = (solid line), corresponding 



to LSI.O, Aa = 0.187 (dashed line), corresponding to the estimate (123|) of Ref. [23|, and 
Aa = 0.115 (dash-dotted line). As in the case of R{t,s), both LSI.O and LSI.l with Aa = 
0.187 do not provide a good approximation of the actual data already for u — 1 < 0.1 
and 1, respectively, whereas LSI.l with Aa = 0.115 shows visible discrepancies only for 
n — 1 < 4 ■ 10^'^, when compared to the curve with s = 2000, which is expected to display 
smaller corrections to scaling. The result of the integration over time has been to reduce 
by almost one order of magnitude the value of n — 1 below which LSI.l with Aa = 0.115 
does not agree with the numerical data in the scaling regime (from ~ 2 ■ 10^^ for R{t, s) to 
^ 2 • 10^^ for X/ (^; s)). By looking solely at Fig. [3] one would naturally conclude that the 
scaling behavior of X/ (^) s) is actually described by LSI.l with Aa = 0.115: The corrections 
for M — 1 < 4 ■ 10^^ seemingly vanish upon increasing s. 

C. Integrated global response functions 

We finally consider the global intermediate response function Xi i'^^s) (see Eqs. (1141) ) 



which was introduced and studied in Ref. 



29| in order to highlight discrepancies between the 



numerical data and the prediction of LSI.O in the case of the the two- and three-dimensional 



Ising model with Glauber dynamics. The same data presented in Ref. 29|] have later on been 



found in agreement with the prediction of LSI.l 23|]. In terms of the quantity ( TT9|) measured 
in the simulations, Xi {^i s) is given by Xi v-i s) = Xk=oit, [s/2, s]) and it is therefore a space 
integrated, or global, quantity. In fact, Xlr -^»'(^'^') ^^ ^^^ response function of the total 
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magnetization. 



(g)/ 



In this Subsection we focus on the comparison between X/ (^) ■^) ^^id LSI, bein 



Ef the one 



29|. 



with the corresponding field-theoretical predictions [3J] already discussed in Ref. 

According to general scaling arguments (see, e.g., Ref. |3|, |29|), the global response func- 
tion Rk^oit, s) is expected to scale as 



Rk=o{t,s) = Ars~'~''+''/^u\u - ly'-'^+'^/^J^niu) 



(39) 



where the non-universal constant Ar is fixed by the condition J-'r{u ^ 1) = 1, and the 
function J^r{u) is universal. Within the RG approach one additionally finds that Tr{\-^) ^ 
and finite., in analogy to Eq. ([7]). As a result of the factorization flT3l) of the space dependence 
of Rr{t,s), the scaling behavior of the Fourier transform Rk=o{t,s) = j drRr(t,s) is given 
by Eq. ([39D with 

(40) 



J^^^^'-'\u)=u-^'^iu 



\Aa 



i.e., it has the same expression as in the case of R{t,s) but with the formal substitutions 
Ar I— > Ar = y4ij Jdx$(x) and a \-^ a — d/z. Note that, as in the case of the local re- 
sponse function, LSI.l and scaling arguments predicts two qualitatively different behaviors 
of !Fr{u -^ 1+). The scaling form of Xi {^^ ■^) is easily calculated on the basis of Eq. 



yielding 



where 



xf^(t,s) = ^«t-'^+'^/V-Vx,g 



.u 



fxM 



dxx 



1/2 



2; \ —1 — a+d/z 
1 Tr{u/x), 

uJ 



(41) 



(42) 



such that /x,g(w ^ 1) = Kg. In particular the prediction of LSI.l is readily found as (see 
Eqs. ([33]) and ^ 



f- 



(LSI.l), 



U] 



dxx 



1/2 



nQ\ —l—a+d/z+Aa 



U 



iN-a /(LSI.l), 



7- 



X,l 



U 



a\-^a—d/z 



(43) 



As in the previous cases, in order to compare the numerical data with the prediction of LSI.l 
we consider the quantity 

g.^,,{t,s)^t'^-"-^u^-\f\t,s) (44) 

which, according to Eq. fHTj) is a function oiu = t/s: 



9xAt^^) =-^RfxM 



(45) 
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such that A^^g = \im.f^sg^^g(t, s) = ArKq. Differently from the previous cases, the behavior 
of g^^g(t,s) for t ^ s is predicted to be qualitatively the same both by scahng arguments 
and LSI.l: The integral in Eq. (H2]) is finite for u —^ 1, whatever the actual behavior of 
J-'r{x/u -^ l"*") is, as long as J^r{u -^ 1+) ~ (n — 1)°" with a > a — d/z, which is the case 
here. In Fig. H] we ploto^g(t, s) for different values of s, which are comparable to those 



considered in Ref. 123, 



291] but generally smaller than the values which Figs. |2] and [3] refer 



to. The data sets corresponding to s = 200 and 400 collapse to a reasonable extent, whereas 
0.24 , , 
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FIG. 4: Scaling function g^^g{t,s) (see Eq. (HH)) as a function of u — 1 (u = t/s > 1) for s = 50, 
200, and 400. The solid line is the prediction of LSI.O whereas the dashed and dash-dotted lines 
are those of LSI.l with Aa = 0.187 and Aa = 0.115, respectively. The two additional dashed 
lines represent the upper and lower bound to Aa provided in Ref. 23(] and resulting from the fit 



of Xi \ti s) with the prediction of LSI.l. In all these cases 6 = 0.38, a = 0.115, and Ar = 0.3895 
have been used. 



the set corresponding to s = 50 still departs from the previous two because of corrections to 
the scaling behavior. For larger values of m — 1, instead, statistical errors increase due to a 
poor statistics. However, our point here is not to determine the scaling function with high 
accuracy, but compare its features to those of the prediction of LSI.l. As it was already 
clear from the numerical results presented in Fig. 2(a) of Ref. 29|| (which Fig. |4]has to be 
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compared to [45|), gx,g{t-, s) does not vary too much as a function of u, being almost flat for 
0-3 ^ ^ ~ 1 ^ 10 while slightly increasing upon decreasing n — 1 < 0.3 (about 2% in the 
interval 0.05 ^ m — 1 ^ 0.3) and decreasing upon increasing m — 1 > 10 (by less than 2%). 
Analogous behavior has been observed in the three-dimensional Ising model (see Fig. 2(b) 
in Ref. 29||). From the large-M behavior of gx,giti s) we estimate the non-universal constant 

^^,g = 0.221(1). (46) 

(Which has actually been obtained on the basis of the data with s = 200 and u > 3.) In Fig.H] 
we also plot ARfx,g ' (u) where fx,g ' (u) in given by Eq. (H2|) and Ar = Ax,g/K,e has been 
fixed to 0.3895. In particular we report the predictions for Aa = (solid line), Aa = 0.115 
(dash-dotted line) and Aa = 0.187(20) (dashed and dotted lines). The prediction of LSI.l 
with Aa = 0.187 provides a quite goo d fit of the actual scaling behavior, in complete 



agreement with what found in Ref. 231], where this value of Aa was determined as the 
best-fit parameter. By looking at the integral in the expression of fx,g ' (u) (see Eq. (1431) ) 
it becomes clear why fx,g ' (u) is so accurate in reproducing the almost constant behavior 
of the numerical data for g^g(u) as a function of u: For Aa = Aao = 1 + a — d/z, fx,g ' (u) 
becomes constant. In two dimensions one finds Aao = 0.192(2) which is within the interval 
determined in Ref. 



23l |. On the other hand, both the prediction of LSI.l with Aa = 0.115 
— which was providing the best fit to R{t, s) — and of LSI.O are quite far from the actual 
numerical data. In particular, the fact that LSI.l with Aa = 0.115 well describes both R{t, s) 
and Xi (t^ s) but not X/ (^7 s) suggests that the factorization of the space dependence (fT3|) 
might not be a good approximation of the actual scaling behavior. 

On the basis of Figs. [2], [3], and H] we can conclude that possible integrations of local 
quantities over time and over space strongly affect the agreement between numerical data and 
the corresponding prediction of LSI.l. Accordingly, the exponent Aa which characterizes 
such a prediction can not be regarded as an additional model-dependent exponents but 
simply as a fitting parameter the optimal value of which depends on the specific quantity 
considered. Indeed, LSI.l with Aa = 0.187(20) describes correctly the behavior of X/ (^; s) 
but it fails in the case of the related R{t, s). Conversely, LSI.l with Aa = 0.115 describes 
very well, for t/s — 1 > 10~^, the behavior of both R{t,s) and x}- (t, s) but it fails in the 
case of the global response function X/ (^5 s)- 
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IV. QUENCHES TO Tf < T^ 

In Sec. [T]we have pointed out that the prediction (fT6|) of LSI.l with Aa 7^ does not 
become stationary in the STSR, in contrast to what the response function actually does. 
Even arguing that the validity of LSI.l might be restricted only to the aging regime does not 
help for quenches at the critical point: The multiplicative structure in Eqs. ([8]) and (TTOj) — a 
direct consequence of scale invariance — leaves no room for LSI.l to describe R^'^^'> without 
spoiling the quasi-equilibrium behavior of R{t, s) in the STSR. For Tf < Tc, however, the 
crossover between i?'^'^'^^(r; Tf) and R^^^\t,s) does not occur multiplicatively but additively, 
as we shall discuss shortly. Then, the criticism moved to LSI.l for critical quenches does 
not extend to below T^. 

Indeed, in quenches below Tc, two-time quantities such as the autoresponse function 
naturally split as [l|, |46|, |47| 

R{t,s) = i?("'i) (r ; Tf ) + R^^^^ (t, s) . (47) 



This additive structure is a direct consequence of the phenomenology of phase-ordering 19| , 
where a patchwork of compact growing domains, internally in equilibrium, evolve through 
defect displacement and/or annihilation. Fast equilibrium fluctuations inside domains are 
responsible for R^'^'^\t; Tf), whereas the slow degrees of freedom associated to the dynamics 
of defects give rise to R^^^\t,s). One easily verifies that Eq. (1471) renders Eq. ([3]) in the 
STSR when R^^^^ is given by Eqs. ([5]) and flTTl) (i.e., by LSI.l): In this sector, in fact, 
u = 1 + t/s -^ 1 and therefore R^^^\t,s) ~ g~Aa^~i-a+Aa -^j^j^}^^ jf /\q, > q^ becomes 
negligible for s — > 00 in comparison to R^°'^\t] Tf). In the aging regime, instead, the opposite 
is true: i?'^'^'^^(r; Tf) for m > 1 becomes negligible compared to R^^^\t,s). In scalar systems, 
which we are interested in here, this is due to the fast (exponential) decay of R^'^'^^r; Tf) over 
a timescale Teq{T{) ~ ,^^»(Tf), where '^(Tf) is the correlation length of the equilibrium state 
at T = Tf with broken symmetry. In the aging regime one eventually has r ^ Teg(Tf) and 
hence R^'^'^\t;T{) becomes (exponentially) small compared to R^^^\t,s) for n > 1, yielding 
R{t, s) = R^^^\t, s), i.e., the response function is very well approximated by R^^^\t, s). 
Now, as stated in the Introduction, in the analytically solved instances of phase ordering 



such as the d = 1 Ising model 10, llll] quenched to Tf = [12] and the time- dependent 



Ginzburg-Landau model within the Gaussian auxiliary field approximation 
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15|, 



Eq. (ITTl) holds |48|. This result has also been verified with high accuracy via numerical 
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Bis 



where 



simulation of the Ising model in 2, 3 and 4 dimensions with Glauber dynamics 
Zg = 2, and of the one-dimensional Ising model with Kawasaki dynamics [37|, where Zg = 3. 
Therefore, in the case of subcritical quenches the form of R^^^^ predicted by LSI.l would fit 
in the analytical and numerical results with 

Aa = 1/zg. (48) 

It must be stressed that contrary to the case of critical quenches, Aa 7^ is not in conflict 
with the time translational invariant behavior of the response function at short times because 



of the additive structure of Eq. (H71) . Note, however, that numerical results for x(t, [C 
and Xk=o{t, [0, s]) in the two- and three-dimensional Ising model with Glauber dynamics 



s 



quenched to Tf < Tc, are seemingly consistent with Aa = (see Ref. 27| for more details). 



llTl, 



This conclusion might be biased by the fact that these quantities, at variance with R{t, s) 



studied in Ref. 



are affected by quite large corrections to scaling 



0. 



V. SUMMARY AND CONCLUSIONS 



We have studied in detail, via Monte Carlo simulations, the impulse autoresponse function 
R{t, s) [see Eqs. ([T]) and ([2])] and the associated (intermediate) integrated local and global 
response functions Xi (t^^) [Eq. (TT7|) ] and X/ (^5 -s) [Eq- (fT4l) ]. respectively, of the two- 
dimensional Ising model with Glauber dynamics, after a quench from high temperature to 
the critical point. 

We highlighted the universal non-equilibrium scaling behaviors of these quantities and we 
compared them with the predictions derived from the theory of Local Scale Invariance in its 
two recently proposed versions, referred to as LSI.O [Eqs. (fT2|) and (fT3|) ] and LSI.l [Eqs. (TT5|) 
and flT3l) ]. The scaling form of R{t, s) predicted by LSI.l depends on the additional parameter 
Aa, compared to the prediction of LSI.O, which is actually recovered for Aa = (see Sec. [T]). 
On the basis of our numerical analysis we conclude that: 

(i) TcR{t, s) becomes — as generally expected — time-translational invariant in the short- 
time separation regime which is accessed by increasing s ^ Tmicr while keeping fixed 
T = t — s (see Fig. [1]). This behavior is consistent with the prediction of LSI.O but not 
with the one of LSI.l (with Aa ^ 0, see Eq. (fTSl)). 
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(ii) The scaling behavior of R(t, s) (see Fig. [2] and Eq. (^^) is correctly captured neither 
by LSI.O nor by LSI.l. However, LSI.l with Aa = 0.115 yields a good fit of the actual 
data only for t/s — 1 > 2 ■ 10^^. Figure[2] — properly normalized by Aji (see Eq. (1281) ) 
— provides a very accurate numerical determination of the universal scaling function 
of R{t, s) in a wider range of times compared to previous numerical studies. 

(iii) Apparently, the scaling behavior of X/ (^i -s) (see Fig. [3] and Eq. ( l35l) ) is correctly 
captured by LSI.l with Aa = 0.115. The small discrepancies which are visible for 
t/s — 1 < 4 ■ 10"'^ can be due either to non-universal correction to scaling (in which 
case they should disappear upon increasing s) or to the fact that the actual scaling 
function of X/ (^; s) deviates from LSI.l at smaller values of t/s — 1 (as in the case of 
R{t, s)). The data at our disposal do not allow one to discriminate between these two 
options. 



(iv) In agreement with previous studies 29|], we find that the scaling behavior of X/ (^) s) 
(see Fig. H] and Eq. ( l44l) ) is encoded in a scaling function g^^g(t,s) which varies less 
than 4% in a wide range of times. Although the quality of the data presented in 
Fig. m does not allow an accurate determination of the scaling function, it is sufficient 
in order to conclude that both LSI.O and LSI.l with Aa = 0.115 have no chance to 



describe it correctly, whereas — as noted in Ref. 23| — LSI.l with Aa = 0.187(20) 
actually provides a better approximation. 

The evidences summarized above lead to the conclusion that the scaling forms predicted 
analytically by LSI.l for the various quantities related to Rr{t,s) provide, in some cases, 
accurate fitting forms but cannot be considered to be exact, at least for quenches at the crit- 
ical point. It would be interesting to compare our numerical determination of the universal 
scaling function of R{t, s) (see Fig. [2]) with the predictions of different analytical approaches 
such as, e.g., field-theoretical ones. 
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